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Abstract 



(N 

We propose an improved parametric form for the equation of state of three- 
l/-} \ dimensional O(N) spin systems. The proposed form is a series expansion with 

two sets of terms, which contribute (mainly) separately to the description of the 
high- and low-temperature regions of the phase diagram. Our goal is a better 



^^ \ description of the low-temperature phase at zero magnetic field (i.e. the coexis- 

tence line), characterized by singularities induced by Goldstone modes. We test 
our proposed form by comparison with existing Monte Carlo data for the N = 4 
case, which is of interest in studies of the QCD phase transition and for which 
the Goldstone- mode effects are quite pronounced. We find that the description of 
the numerical equation of state is indeed improved with respect to other fitting 
forms. In all cases considered we determine the coefficients nonperturbatively, 
from fits to the data. As a consequence, we are able to obtain a very precise 
characterization of the pseudo-critical line for the model. 



1 Introduction 

The O(N) (or, more specifically, the N- vector) spin models correspond to a general- 
ization of the Ising model to the case of the continuous symmetry of rotation. The 
spin variables Sj are taken as vectors on a sphere of unit radius in an iV-dimensional 
space. We consider N > 2. The Hamiltonian is denned in terms of the scalar product 
of nearest-neighbor spins on a three-dimensional square lattice as 

pn = -j J2 Si-Sj - h ]rs,, (i) 

where J > represents the ferromagnetic coupling and H the external magnetic field. 
These models are of general interest for the statistical mechanics of phase transitions 
PP. The N = 2 case (also known as the XY model) describes the superfluid transition 



in liquid helium and the N = 3 case corresponds to the classical version of the Heisen- 
berg model for ferromagnets. 1 Moreover, it is believed that the N = 4 case describes 
the chiral phase transition in finite-temperature QCD with two degenerate light quark 
flavors, which makes this class of models interesting to high-energy physics as well. 2 
In this case, the magnetization and the magnetic field of the spin model correspond 
respectively to the chiral condensate and to the quark mass for the QCD analogue of 
the transition 0E1I11- 

The O(N) symmetry is exact in the Hamiltonian for H = 0, just like the reflection 
symmetry for the Ising model. The main difference with respect to the Ising case is the 
possibility of configurations where the spins are locally aligned but for long distances 
this alignment is lost, yielding a null average for the magnetization. Such configurations 
- called spin waves — possess arbitrarily low energy and tend to destroy the order of 
the system even at low temperatures. In fact, as opposed to the Ising model, the O(N) 
models do not display a phase transition with spontaneous magnetization 3 in d — 2. 
In d — 3 a phase transition occurs, with the presence of spontaneous magnetization 
below the critical temperature. The breaking of the (continuous) rotational symmetry 
at low temperatures, signaled by the spontaneous magnetization, is associated with 
Goldstone modes, the spin waves. These modes cause the divergence of the zero-field 
susceptibility not only at the critical temperature, but for the entire low-temperature 
phase [HE]- Note that the magnetic field defines a privileged direction in spin space and 
the magnetization M is the expectation value of the spin component along H. There 
are thus N — 1 massless Goldstone modes, corresponding to the iV — 1 transverse spin 
components. 

Spin models in the O(N) class have been extensively studied using analytic and 
numerical methods (see [U] for a recent review). In particular, the nonperturbative study 
by Monte Carlo simulations is very efficient for these models due to the Swendsen-Wang 
cluster algorithms 7j, which can be applied to the continuous-spin case by means of 
the embedding technique introduced by Wolff jH]. This study is important to test the 
perturbative predictions and to investigate cases for which these predictions are not 
available, or cannot be done with great accuracy. These problems include properties 
of the models in the presence of magnetic field and the direct calculation of long- 
distance observables such as the correlation length. For example, the predicted singular 
behavior of the longitudinal susceptibility for vanishing H - - mentioned above and 
induced by Goldstone modes at low temperatures -- was directly observed in Monte 
Carlo simulations of the cases N = 2, 4, 6 respectively in Refs. 0, [HI] and [TT] . 

Here we consider the determination of the magnetic equation of state, which gives 
the relation between applied field, temperature and magnetization for the system. The 



lr The N — and N — 1 (the Ising model) cases, not considered here, correspond respectively to 
models for the statistical properties of long polymers and for the liquid- vapor transition in several fluid 
systems. 

2 Two-dimensional O(N) models are also of indirect interest in quantum field theories, as toy models 
for asymptotically-free gauge theories. 

3 For the case N = 2 there is a phase transition of the Kosterlitz-Thouless type, without spontaneous 
magnetization. 



equation of state has been determined perturbatively for general N by e-expansions (see 
PU Chapter 29] and references therein) and for the cases N = 2, 3, 4 by matching a 
high-temperature expansion (with coefficients obtained from perturbation theory) to a 
parametric form incorporating the leading Goldstone-mode behavior jHJEl- Of course 
it is interesting to compare these expressions to Monte Carlo data for the equation 
of state. One can also test the various forms used in the perturbative expansions 
(or new proposed forms) by fitting them to the Monte Carlo results and obtaining 
nonperturbative coefficients. This has been done (see e.g. [10J) using an interpolation 
of the low-temperature (Goldstone-mode) form derived in ^3] with a high-temperature 
form determined by analyticity conditions. This method has the advantage of a clear 
low-temperature form, with several orders in the Goldstone-mode expansion, but has 
the disadvantage of needing an interpolation with the high-temperature form. 

In the present paper we carry out fits using instead a variant of Josephson's parame- 
trization PIE]; a polynomial parametric representation for the equation of state. The 
resulting representation is valid above and below the critical temperature and automat- 
ically satisfies the analyticity conditions mentioned above. In addition to the leading 
(multiplicative) Goldstone-mode contribution, we consider explicitly the higher-order 
terms, which are important in the low-temperature region. Our proposed form con- 
tains two sets of coefficients, which will be separately more relevant for the description 
of the high- or low-temperature regimes. We argue that the use of this double set of 
coefficients enables a better characterization of the two regimes, leading to better fits 
in the comparison with numerical data. This claim is verified by an application to 
the N = 4 case, for which the Goldstone-mode effects are fairly high, using the data 
reported in JU|. As mentioned above, this case is of interest for comparison with data 
from numerical simulations of the phase transition in two-flavor QCD. In particular, the 
prediction of universal behavior in the 0(4) class has been confirmed for lattice-QCD 
data in the Wilson- fermion case |15j . but not for the staggered- fermion formulation, 
which is believed to be the appropriate formulation for studies of the chiral region. (At 
the same time, some recent numerical studies suggest that the transition may be of 
first order 16 .) We plan to extend our analysis to the N = 2 case, for which we are 
generating new data [T7j . 

The paper is organized as follows. In Section El we describe the usual parametric 
representation for the equation of state, as well as our proposed form. In Sections El and 
0]we consider the determination of important universal properties that can be obtained 
from the equation of state: some critical amplitude ratios and the characterization of 
the pseudo-critical line (respectively in Section 01 and in Section HJ). Finally, in Sections 
03 and El we present our results and conclusions. 

2 Scaling equation of state 

The magnetic scaling equation of state is given jTJ Chapter 29] by 

h = M s f{t/M 1/t3 ), (2) 



where t and h are the reduced temperature t = (T — T c )/T and magnetic field 
h = H/H . We fix the normalization constants T and H by requiring unit critical 
amplitudes in the behavior of the magnetization along the coexistence line (given by 
t — > 0_, h = 0) and along the critical isotherm (given by h — > 0, t = 0), corresponding 
respectively to M = {-tY and M = h 1 / 5 . 

The equation of state can also be written as 

y = fix) , (3) 

where 

y = h/M s , x = tjM ll<i . (4) 

Note that the coexistence line and the critical isotherm are given respectively by x = — 1 
and x = 0. The corresponding normalization conditions are thus 

/(0) = 1, /(-l) =0. (5) 

For large values of x (i.e. in the high-temperature region of the phase diagram) the 
behavior of f(x) is described by Griffiths's analyticity condition P 

oo 

/(ar) = $> n ^- 2 ("-^ . (6) 

71=1 

As said in the Introduction, at low temperatures there appear divergences in the 
zero-field magnetic susceptibility due to transverse fluctuations from the massless Gold- 
stone modes [TJ El CHI UH1 120] • To leading order the divergence in the longitudinal sus- 
ceptibility is proportional to h^ 1 ^ 2 and the equation of state has the leading behavior 

f(x) = y oc (1+a;) 2 (7) 

for x — > — 1. We note that the Goldstone-mode divergences cancel out and the equation 
of state is divergence free. This is seen order by order in the e-expansion [21] and 
for fixed- dimension perturbation theory 22. This is also observed nonperturbatively 
in the Monte Carlo data (see e.g. ^0])- The corrections to the leading behavior are 
incorporated explicitly in the expression proposed by Wallace and Zia [T3] . which is 
inferred from the e-expansion for the equation of state deduced in [21] . For d = 3 the 
expression corresponds to an expansion in powers of y 1 ^ 2 

xi(y) + 1 = (ci + d 3 ) y + c 2 y 1/2 + d 2 y 3/2 + ■■■ (8) 

This form describes well the Monte Carlo data from the low-temperature region until 
around the critical temperature. The coefficient associated with the H~ 1 ^ 2 divergence 
of the susceptibility for H — > is Z 2 - Note that the expression of c 2 derived in JH] 
increases with N, i.e. models with larger iV should display stronger Goldstone-mode 
effects. 

In References jD] , JO] and ^1] the Monte Carlo data for the equation of state have 
been fitted to the expression 

x(y) = x x {y)^j-^ + x 2 (y)^— , (9) 

yd + y y$ + y 



M 


= m R^m(9) 


t = 


r(i- e 2 ) 


H 


= h R fi5 h(9) . 



where x\(y) is given in Eq. (jHJ) above and 

x 2 (y) = ay 1 '"' + by {1 ~^ h (10) 

corresponds to the first two terms of Eq. (JHJ. (The parameters y$ and n are chosen to 
ensure a good interpolation.) This interpolation of low- and high-temperature behaviors 
describes well the data, but of course it would be nicer to have a form valid in both 
regions, such as the parametric form introduced in j!4j . This type of form is described 
in the next section. We then comment on the previous use of this parametrization and 
propose a new variant that is especially well suited for fits. 

2.1 Parametric representation 

Let us consider the polynomial parametric representation introduced in [Hj, in which 
one writes M, t and H in terms of the variables R and 9 (see e.g. [D E]) 

(11) 
(12) 

(13) 

Here m(9) and h(9) are odd functions 4 of 9, regular at 9 = and 9 = 1. This ensures 
that Griffiths's analyticity conditions are satisfied. The coexistence line is given by 9q, 
the smallest positive zero of h(9). [From Eq. (|12j) it is clear that 9$ must be greater 
than 1.] Without loss of generality, we may take m{9) = 9. The equation of state then 
becomes 

x = W^{i) • (14) 

m = r-M . d5) 

The relation between x and 9 is shown schematically in Fig. ^ together with the re- 
spective ranges considered. Note that we must have 9^ < 1/(1 — 2/3) for the above 
mapping to be invertible [21]. For the 0(4) case (3 ~ 0.4 and we have roughly 9$ < 5. 

With the parametrization (fTT |) -(fT3* )) . the singular part of the free energy T s can be 
written as 

T s = h m R 2 ~ a g(9), (16) 

where the function g(9) is the solution of the first-order differential equation 

(1 - 9 2 )g'{9) + 2(2 - a)9g{9) = \{1 - 9 2 )m'{9) + 2(39m{9)} h{9) (17) 



4 The function h{9) should not be confused with h, the normalized magnetic field introduced in 
Section [3 
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Figure 1: Schematic representation of the relation between x and 0. From left to 
right, the ticks correspond respectively to the coexistence line, the critical point and 
the high-temperature (zero- field) limit. 



that is regular at = 1 [2B]- (This follows from the relation H = dT s jdM , where the 
derivative is taken at fixed t.) 

We note that the parametrization above was first used in perturbative studies of 
the equation of state for the Ising model 24J. We discuss below its application to the 
iV-vector (Goldstone-mode) case. 

In accordance with Eq. (J7J) the leading behavior for 9 — ► 9q must be 

h{9) -»■ (0 O - 9) 2 for 9 -> O ■ (18) 

This combined with the requirement that h(9) be an expansion in odd powers of 9 
suggests the general form 

h{9) = 9(l-9 2 /9 2 ) 2 (l + jr Cl 9A . (19) 

This form is used in jHl^j, in their "scheme B". They also define another scheme with 
a similar expression for m(9). In both cases the differential equation becomes 

3+n 

(1 - 9 2 )g\9) + 2(2 - a)9g{9) = 9 £ a { 9 2i , (20) 

i=0 

with coefficients c^ depending on the exponent (3, the root O , the coefficients q and 
on the scheme considered. One can easily check that the solution of the differential 
equation that is regular at = 1 is given by 

3+n i ^ £ Q2{i-k) /j _ q2\* 

9{9) = " S S i(i-k)l( a -2)...( a -2 + k) ■ (21) 

Clearly, this solution is a function of the values of 0o and of the parameters Cj, i = 
1, . . . , n. In [HI E] the authors have considered the cases to = 1,2, with parameters 0o, 
Cj obtained from perturbation theory. We comment on their results in Section 15.21 

In the next subsection we introduce a more specific parametric expression for h(9), 
as a combined expansion around = and around = O , in order to isolate the two 
regions of the phase diagram. 
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2.2 Improved parametric form 



We consider here a variant of the parametric function h{9) above 



h(9) 



9 1 






i + E c ^ 



2; 



1 + H 1 




1 



2\ J 



/ 1 + E^i 



(22) 



For consistency, the expression is normalized so that the contribution from the d/s is 
equal to 1 at 9 = 0. This normalization factor does not affect the equation of state, 
since h(6) enters in f(x) only as a ratio. Let us note that this is still an odd function 
of 9 and is equivalent (as an expansion in 9) to equation (|19|) with a rearrangement of 
terms. In particular, we may compare the series using only terms with q coefficients 
(i.e. with all dj = 0) with the one using only dj coefficients. The relation between the 
two cases is given by 



C; 



III 

-d*eHH 



9f (l + Y,*i 



(23) 



Thus, a single coefficient q corresponds to a sum of d/s. Conversely, if we considered 
an expansion around 9 ~ 9q each dj would correspond to a sum of q's. Since the roles 
of the terms of the two series are different, considering both series may be important 
for truncated sums, such as the ones we use for the fits. 

The consideration of two types of coefficients (q and dj) is done for gaining better 
control over the description of the two distinct regions of the phase space, the low- 
and high-x regions. In fact, although the two sets of coefficients give rise to a similar 
expansion in powers of 9, the determination of the single coefficients (by fits to the 
numerical data) is more stable when each of the two regions is separately associated 
with a set of coefficients. More precisely, since we write the series as a product of two 
sums of terms corresponding respectively to an expansion around the high-x region 
(9 ~ 0, coefficients Cj) and the low-x region (9 « 9q, coefficients dj), we can expect fits 
of the data for each of these two regions to be more sensitive to the corresponding set 
of coefficients, since the other set's main contribution will be a constant. An indication 
of this property can be seen from a "quick" expansion of the parametric form in powers 
of e ~ for the two regions 



h{e) 



e 1- 



2_£2 

9 2 



1 + cie' 




9 2 



/(1 + E^I (24) 
3=1 



«*-.>- ¥K 



1 + E^M 1 



4 = 1 



2ie 



(1 + o?! 2 e/6"i 



1 + J2dj 

3=1 



(25) 



Note that we show only the leading order from each multiplicative contribution. It is 
interesting that the correction to the leading behavior is C(e 2 ) in the first case and 0(e) 
in the second. 
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To be more precise, we can associate the various q, dj coefficients with the coeffi- 
cients in the separate expressions for the low- and high-x regimes used in [TO], respec- 
tively equations (jHJ) and (fTUj) above. For the high-a: behavior we expand f(x) around 
9 — * (corresponding to large x). We start by writing 9 as a function of x and inverting 
Eq. (|14|) consistently in powers of 9 (correspondingly in powers of x' 1 ^). We then get 



where 



m{6) 


= 9 = (Ax)- 13 ' 


h(9) 


= (Ax)- [l + B 




A = (91 - 1) *o 



13 (Ax) 



-2/3 



0(x- 4f3 ) 



-2/3 



+ 0(x 



-4/3 > 



w 



B = Cl -- 2 -(3-- Ejidj 

9 Q 9 Q 1 + J2j dj 



The equation of state then becomes 

(A X y 



u 



h(i) 



1 + ((35 + B) (Ax)- 2 " + 0(x 



-4/3 \ 



(26) 

(27) 

(28) 
(29) 

(30) 



This expression is of the form (JHJ) and can be inverted and compared to Eq. (jlUj) , giving 



j Hi)] 

_l-20 



1/0 



7 



((35 + B) 



(31) 



(32) 



Note that the leading coefficient a contains 9q and sums of the coefficients q and dj, 
whereas the expression for the next orders will contain isolated contributions from the 
Cj's (e.g. the coefficient c\ in the expression for b) but not from the dj's, which appear 
always as a sum. 

Analogously, for the low-x region we expand the expressions of 9, h(9) around 
9 — * 9 (corresponding to x — » — 1) and substitute the results into the expression for 
f(x). Defining 

9 = 9 (l- e) (33) 

we write x as a function of e, invert this expression to get e(x) and then obtain h(9) in 
terms of x, as done above for the large- x case. The expressions are 



h(9) 



(l+x\ 




{ A' ) 




n( l + 


X 



B' 



1 + x 
A' 



IB 



a 



C_ 
I 7 



.4' 



E 



+ F 



2B'^ 

3B'E^ 
A> 



A' 2 

l + x 
A' 

2C 5B' 2 

+ ~A^ 



l + x 
A' 



+ 



(34) 



l + x 
A 1 



+ 



(35) 



with 



A' = 
B' = 
C = 



_26l_ 
9 2 -1 

01 



ft 2 



1 



-<MH8 +2 , 



ol 



o 2 



1 
1/F 



D = 4^ (l+J^Cie^/ll + ^dA 



E 



2di - 2 



1 + E,Q^ 



F S 4^-5d 1 + g + 4 ^-, dl I £ ^ + " ;(2; - 1) ^" 



2; 



1 + E,Q^ 



l + EiCi^ 



(36) 

(37) 
(38) 

(39) 
(40) 
(41) 



The equation of state then becomes 



D6 5 f\ + x 
h(l) V A' 



l+ ( E -™ + 6) n ' 



A' 



A 1 



G 



1 + x 
A' 



+ 



where 



G 



5(5+1) r ^ „ 355' + 35'£ + 2C 55' 2 
+ 5 h -- h -g- . 



(42) 



(43) 



2 A' A" 

This form may be inverted to give an expression of x as a series of powers of y 1 ^ 2 as in 
Eq. ©. We obtain the coefficients 



c 2 = A' 



D 



1/2 



?2 
c 2 



^+^3 = ~^[E 



2A> 



2B' 



+ 5 



do 



2A! 



Ue- 2 4 + s 

4 I A' 



G 



(44) 
(45) 

(46) 



We see that in this case it is the dj coefficients that appear as single contributions, 
whereas the q's appear always as sums. 

Thus, the qualitative feature observed in Eqs. (J2*lj) and (|2*3|l is confirmed by a more 
careful expansion, i.e. the c's are more relevant for the high-x region and vice-versa. 
This will also be seen directly from the fits in Section 15.11 All calculations above were 
checked using Mathematica. 
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3 Amplitude ratios 

Just like other critical properties of statistical systems (e.g. critical exponents), certain 
ratios of critical amplitudes are universal [25]. The amplitude ratios are taken as dimen- 
sionless combinations of critical amplitudes above and below T c for various quantities. 
For example, for the singular part of the specific heat one has 

C H = A ± \t\' a , t^±0, (47) 

where t oc (T — T c ). The ratio A + jA~ is then universal. Similarly, by considering the 
behaviors of 

• the susceptibility along the critical isochore (t > 0, H — 0) 

X = C + t~^ (48) 

• the magnetization along the critical isotherm (t = 0, H ^ 0) 

M = D- 1/S H 1/S (49) 

• the magnetization on the coexistence line (t < 0, H = 0) 

M = B(-tf (50) 

one may construct the universal ratios 

R c = aA + C + /B 2 , (51) 

R x = C+DcB 6 - 1 . (52) 

These and other universal ratios may be obtained directly from Monte Carlo simulations 
(as done e.g. in J2E]) or indirectly from the equation of state. In the case of the Josephson 
parametrization discussed above, the universal amplitude ratios of quantities defined at 
zero momentum are given in terms of g(8) by [21] 

A+/A- = {91 - If- M- , (53) 

R c = -«(l-«)(2- a ) (gg ~,^ (Q) , (54) 



(# 2 -l)^'(0) 



R x z : 77?2 TvTTTTTa • ( 55 ) 



Note that in order to evaluate the function g(8) one has to solve the differential equation 
(JT7J), i.e. determine the solution (|2*T|). 

Our results for the above ratios are reported in Section 15.31 
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4 The pseudo-critical line 

Another important property that can be extracted from the equation of state is the 
characterization of the so-called pseudo-critical line, defined by the points where the 
susceptibility \ shows a (finite) peak for H ^ 0. This corresponds to the rounding of 
the divergence observed at the critical point, i.e. for H = and T = T c . More precisely, 
one looks for a peak in the scaling function of the susceptibility, given by 



M 



h 1/5 f M (z) =» X 



dM 
~dH 



h 1 ' 6 



-i 



#o 



/xW 



where 



z = t/h 1 ^ . 
Clearly, at each fixed h the peak in x is given by t p = Zph 1 ' 136 , and we have 



M n 



h^fuizp), H oXp = h 1 ' 5 - 1 f x (z p ) . 



(56) 
(57) 
(58) 



Thus, the behavior along the pseudo-critical line is determined by the universal con- 
stants z p , fM(zp), fx( z p)- Determining this line is important for systems where a study 
at H = is not possible (and consequently the critical value T c is not known with 
accuracy), such as for the chiral transition of QCD at finite temperature. In fact, the 
knowledge of these universal constants allows an unambiguous normalization of QCD 
data (using the observed scaling along the pseudo-critical line), as done in [27| . 

The pseudo-critical line has been studied for 0(2) and 0(4) models in |2*Hll29j . For 
the N = 4 case, it is found that the susceptibility peaks are given by z p = 1.33(5). Since 
this value is close to the interpolating point of the equation of state in [TJIj, it is very 
important to work with the smooth parametrization considered here, especially when 
using the derivative of /m(^) as i n Eq. (joT?|) below. 

The expression for f x (z) can be easily obtained from the equation of state, given by 
fivi(z) or f(x). Using the original parametrization we obtain 



fx(*) 



1 
6 



fM(z) - -gfM(z) 



£ [/(S)] 



1-1/5 



f35f(x) -xf'(x) 



(59) 



(Note that z = x [f(x)] 1 ^ s '.) In terms of the parametric representation this gives 



fxlfi) 



z{9) 



Mi) 



~ 1/S (2(36 2 + 1 - 6 2 )h(6) 
2f356h{6) + (1 - 6 2 )h'(0) 

-i//m i/e (1 _ g2) 

91 - 1 • 



(60) 
(61) 



Our results for f x (z) and the determination of z p are shown in Section I5~31 
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ol 


C\ 


c 2 


X 2 /d.o.f. 


2.33(3) 






0.50 


2.01(8) 


0.16(6) 




0.43 


1.67(5) 


0.22(6) 


0.18(4) 


0.44 


el 


d\ 


d 2 


x'/d.o.f. 


3.61(4) 


438(1) 




0.59 


7(2) 


-1.4(2) 


-0.075(1) 


0.44 



Table 1: Fits in the high-temperature regime (using x > 0). The values of \ 2 jd.o.j . 
should be taken only as relative measures of the goodness of the fits. The number of 
d.o.f. is 33. 



5 Results 

The fits have been done using a conjugate-gradient minimization [3U| of \ 2 — without 
considering the gradient of the function f{6) — with a numerical inversion of Eq. (|14|) in 
order to find 9 for any given value of x. For the critical exponents we used v = 0.749(2) 
EH and 8 = 4.824(9) j22|, implying the values (3 = 0.386(1), 7 = 1.476(5) and the 
upper bound 6q < 4.38(5). We refer to these values as the first set of exponents. We 
note that the corresponding exponent 8 from (HIj, 4.789(6), has slightly smaller error 
bars. However, we choose to use the one in [32] because it is obtained directly from 
(infinite-volume) simulations at nonzero magnetic field. These two exponents are not in 
agreement within error bars. We will also present below for comparison a few quantities 
obtained using the exponent S from [ST]. We refer to the resulting values as the second 
set of critical exponents. 

The data for the magnetization are taken from [10;. In addition to the statistical 
errors, we have included errors due to the critical exponents, the critical temperature 
and the normalization constants H$ and To. These constants have been rederived using 
the first set of exponents above (with errors), yielding 



H = 4.85(2) , T = 1.055(5) . 



(62) 



The errors reported in parentheses in all the tables below are Monte Carlo (MC) errors, 
obtained with 2000 MC iterations. In particular, in Section Iq~TI we not only vary the y 
variable but also consider the uncertainties in the exponents 7 and 5 appearing in the 
fitting function [i.e. in Eqs. (|T4*j) and ()15[) ]. The same is true for the errors reported in 
Section \5. 21 In Section IB~3*1 the errors comprise the error bars in the input parameters 
and also the errors in the critical exponents. 
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el 


C\ 


c 2 


X 2 /d.o.f. 


1.905(7) 






40.1 


1.09(2) 


-1.18(4) 




27.3 


1.07(1) 


3.3(1) 


-5.1(1) 


28.1 


01 


d\ 


d 2 


X 2 /d.o.f. 


3.85(4) 


-4.0(2) 




20.2 


2.69(2) 


154(4) 


-111(2) 


18.7 



Table 2: Fits in the low-temperature regime (using x < 0). The values of x 2 /d.o.f. 
should be taken only as relative measures of the goodness of the fits. The number of 
d.o.f. is 37. 



01 


C\ 


c 2 


c 3 


c 4 


X 2 /d.o.f. 


1.955(7) 










31.5 


1.614(7) 


0.58(3) 








19.6 


1.392(5) 


-0.06(1) 


0.80(3) 






18.1 


1.247(6) 


1.6(2) 


-2.8(3) 


2.7(2) 




17.6 


1.170(3) 


-0.7(1) 


6.8(4) 


-11.6(8) 


7.2(5) 


17.4 



Table 3: Fits using only q terms and the whole set of data. The values of x 2 /d.o.f . 
should be taken only as relative measures of the goodness of the fits. The number of 
d.o.f. is 69. 



5.1 Fits 

As a first step, we tried to fit the data separately in the high- and low-x regimes using 
only Q or only dj parameters, in order to confirm that the q's are more important at 
high x and the d/s at low x, as suggested in Section l2~2l As one can see from Tables 
H] and 121 this is indeed the case. At high x the fits using q parameters work better 
that the fits using dj parameters, as can be seen in the case with one parameter plus 
0q. When using two parameters plus 6q the values of x 2 /d.o.f . obtained in the two 
cases coincide, but in the case with d\ and d 2 one obtains the unphysical value 9^ « 7. 
Moreover, if one tries to do a fit using 6q, c\, c 2 and d\ the fit is not better than the 
one reported in the third row of Tableland the value of d\ is very close to 0. In the 
low-x region the fits using dj parameters work much better than the corresponding fits 
using q parameters. Again, if one tries a fit using 6*q, Ci, d\ and d 2 the result is not 
better than the one reported in the last row of Table 121 Thus, we see clearly that the 
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el 


d x 


d 2 


C?3 


G?4 


X 2 /d.o.f. 


3.99(4) 


-3.5(1) 








12.0 


3.22(4) 


-9.2(9) 


3.9(5) 






10.2 


2.63(2) 


-69(3) 


83(2) 


-36(1) 




10.0 


2.73(2) 


-53(3) 


42(2) 


5.9(2) 


-15.0(9) 


10.2 



Table 4: Fits using only dj terms and the whole set of data. The values of x 2 /d.o.f . 
should be taken only as relative measures of the goodness of the fits. The number of 
d.o.f. is 69. 



fit (ci + dj) 


f) 2 


C\ 


C2 


C3 


di 


x'/d.o.f. 


3+1 


2.16(3) 


0.80(6) 


-0.39(7) 


0.58(4) 


33(6) 


9.8 




el 


C\ 


C2 


di 


d 2 


x'/d.o.f. 


2+2 


2.17(4) 


0.9(1) 


-0.62(7) 


-1.56(4) 


1.15(5) 


9.8 




el 


C\ 


di 


d 2 


d 3 


X z /d.o.f. 


1+3 


2.16(2) 


1.4(1) 


31.2(9) 


-50(1) 


38(2) 


9.8 



Table 5: Fits using 5 parameters and the whole set of data. The number of coefficients 
Ci and dj used in each case is indicated in the first column. The values of x 2 ''/d.o.f. 
should be taken only as relative measures of the goodness of the fits. Here we use the 
first set of critical exponents. (The number of d.o.f. is 69.) 



coefficients q and dj are more relevant respectively at high and low x, as suggested in 

Section O 

As a second step, we checked that the fit of all the data using only the parameters 
Ci does not work very well (see Table EJ) • In particular, even with four parameters Cj 
one cannot get a large improvement in the value of \ 2 /d.o.f., compared to the case 
with only the parameter c\. The situation is slightly better when considering only dj 
parameters (see Table EJ). Notice, however, that we have only a few data points with 
very large x and that the low-x expression used in jTU] describes well the data up to 
x « 2. 

Finally, fits of all the data with both c; and dj parameters (see Table EJ) work 
very well, giving a value of \ 2 /d.o.f . about a factor two smaller than the best result 
obtained in Table El (see last row). In particular, we find it interesting that for the three 
fits considered we get (within errors) the same value for Q\. By averaging over the three 
results we find 

di = 2.16(2) . (63) 

In Fig.|2]we show a plot of the data together with the curve corresponding to the case on 
the second row of Table El We have also tried fits with 6 parameters, without significant 
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Figure 2: Plot of the data together with the fitting curve for the case with coefficients 
Ci, C2 and di, d 2 using the first set of exponents. No errors are shown for the curve. 
Error bars on the data are one standard deviation. 



improvement in the value of x 2 /d.o.f. 

The relatively high values of x 2 /d.o.f. may be related to remaining systematic effects 
in the data. This is especially true in the low-temperature regime, where the finite-size 
effects are very strong due to the effect of Goldstone-mode-induced singularities. It 
would be interesting to test our parametrization using the higher-precision data recently 
produced in [32]. In any case, if we use the second set of exponents above (i.e. with 5 
from J3T]) the values of x 2 /d.o.f. are significantly worse, as can be seen in Table |H1 We 
note that, in order to consider this second set of exponents, we have reevaluated the 
normalization constants Hq and To, the values of x and y and the data errors for this 
case. 



5.2 Comparison with other parametrizations 

We now compare our results with previous expressions for the 0(4) equation of state. 
In Reference [12] , the scheme B considered by the authors corresponds to all dj = 
and only Oq, c\ nonzero. Their values for these coefficients are 



01 



2.4(2) , Cl = 0.065(30) . 



(64) 



Note that their value of 8q is consistent with ours within error bars. Using these 
coefficients as a "fit" of the data (considering the first set of critical exponents above), 
one obtains a x 2 /d.o.f. of 268. We can also use our second set of data to evaluate 
X 2 / d.o.f., but this yields the value 688. 

We also consider the interpolated parametrization introduced in [111], presented in 
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fit (cj + dj ) 


e'l 


Cl 


Cl 


C3 


c 4 


X'/d.o.f. 


4+0 


1.137(2) 


-1.0(2) 


8.5(6) 


-14.6(9) 


9.0(5) 


58.1 




el 


Cl 


C2 


C3 


rfi 


xVd.o./. 


3+1 


2.054(7) 


1.17(6) 


-1.05(6) 


1.11(4) 


510(30) 


27.0 




el 


Cl 


C2 


di 


rf 2 


xVd.o./. 


2+2 


2.19(3) 


1.2(1) 


-0.80(6) 


-1.74(3) 


1.34(4) 


26.6 




el 


Cl 


dx 


d 2 


d 3 


X'/d.o.f. 


1+3 


2.22(1) 


0.67(1) 


630(20) 


-920(20) 


560(10) 


29.0 




el 


d x 


d 2 


d 3 


C?4 


xVd.o./. 


0+4 


2.13(1) 


-71(2) 


180(4) 


-193(3) 


75(2) 


26.7 



Table 6: Fits using 5 parameters and the whole set of data. The number of coefficients 
q and dj used in each case is indicated in the first column. The values of x 2 /d.o.f. 
should be taken only as relative measures of the goodness of the fits. Here we use the 
second set of critical exponents. (The number of d.o.f. is 69.) 



Eqs. flEJ), © and (fTU|) above. Using the first set of data, we obtain the high-x coefficients 

a = 1.07(1), b = -0.95(3) (65) 

with x 2 /d.o.f . = 0.52 (cut at x — 1.5). At low x we get 

ci + d 3 = 0.19(1) , c 2 = 0.746(3) , d 2 = 0.061(8) (66) 

with x 2 /d.o.f. = 25.5. Note that the above coefficients are only in partial agreement 
with the values in [TU] and J22], mostly due to the slightly different critical exponents 
considered. We then use these coefficients for the interpolated expression in Eq. 0, 
setting (as in [TU] ) y = 10, n = 3. The resulting 5-parameter fit of the data has 
X 2 1 d.o.f. = 26.2. 



5.3 Universal quantities 

As discussed in Sections E3 and El we use the fits obtained above to evaluate several in- 
teresting universal quantities, such as critical amplitude ratios and the characterization 
of the pseudo-critical line in the phase diagram. 

In Table [7| we show the results obtained for the ratios A + /A~, R c , R x using our 
preferred fits (reported in Table EJ). The three fits give consistent results within error 
bars. Averaging over the three cases yields 



A + /A~ = 1.8(2), 



Rr 



0.26(1) 



R x = 1.10(5) . 



(67) 



These values are in agreement with the ones reported in ^21 Table 3]. (Note, however, 
that our values take into account the errors due to the uncertainty in the critical expo- 
nents.) We also show, in Table El the same quantities using the fits for our second set 
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fit (cj + dj) 


A+/A- 


R c 


R x 


3+1 


1.7(2) 


0.26(2) 


1.11(6) 


2+2 


1.8(5) 


0.26(2) 


1.1(1) 


1+3 


1.8(4) 


0.26(2) 


1.1(1) 



Table 7: Results for the universal amplitude ratios using the fits reported in Table El 



fit (a + dj) 


A+/A- 


R c 


R x 


3+1 


1.6(1) 


0.25(1) 


1.09(5) 


2+2 


1.6(1) 


0.27(2) 


1.1(1) 


1+3 


1.9(3) 


0.22(1) 


1.02(5) 



Table 8: Results for the universal amplitude ratios using fits reported in Table El 



of data (from Table EJ). These results show a little more fluctuation, but are essentially 
in agreement with the ones in Eq. ()67|) above. Note that the ratio R x can also be 
evaluated directly [9 J from the coefficient a in the interpolated form, given in Eq. ([65)1 . 
In this case we get R x = a 
Ref. [32]. 



1.105(15), in agreement with our result above and with 



We now turn to the numerical characterization of the pseudo-critical line (see Section 
EJ). Using Eqs. (|6T)|) and (jfilj) . we draw the parametric plot of the scaling function for 
the susceptibility versus z (see Fig. El). The peak corresponds to the pseudo-critical line 
and can be determined numerically from the two equations by varying 9. The peak 
coordinates thus obtained are reported in Table 01 where we used our preferred fits. 
The values are consistent within errors, yielding 



9 P = 0.587(2) 



1.29(1) 



f x {z p ) = 0.341(1) 



The results are in agreement with previous determinations of z p and /. 



XV^P/' 



(68) 



made in 



Refs. J2H] and [T2|, but our error for z p is much smaller. In Table ITU1 we present these 
quantities in the case of our second set of data. Again, the determinations are in 
agreement. 



6 Conclusions 

We have introduced an improved parametric form for the description of the equation of 
state of 3d O(N) models. This form is based on the parametrization used perturbatively 
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Figure 3: Plot of the scaling function of the susceptibility f x (z) from Eqs. (J5U|) and 
(joTj) , using the fit in the second row of Table El 



fit (a + dj) 


6 P 


Zp 


Jx\ z p) 


3+1 


0.580(2) 


1.33(1) 


0.340(2) 


2+2 


0.589(4) 


1.28(3) 


0.343(3) 


1+3 


0.592(3) 


1.27(1) 


0.339(2) 



Table 9: Results for 9 P , z p and f x (z) using the fits reported in Table El 



in [21] for the Ising model, but takes into account terms associated with the effects 
of Goldstone-mode fluctuations. Such effects are present in O(N) models along the 
coexistence line, i.e. at low temperatures and small magnetic field (or equivalently, at low 
values of the variable x). These new terms are included by means of the dj coefficients, 
associated with an expansion around the coexistence line. (The d/s are considered in 
addition to the usual c, coefficients, related to the high-temperature/high-x behavior.) 
We show that the new parametric form indeed provides a better fit to the numerical 
data as compared to previous parametrizations. In particular, the consideration of the 
dj coefficients is essential for a good description of the Monte Carlo data in the whole 
range of values of x. Also, we were able to verify clearly the different roles played by c, 
and dj parameters in the high- and low-x regions. 

We note that --in the case where all dj = -- our parametrization is equivalent 
to the scheme B discussed in jEJ, used perturbatively by the authors for general 0(N) 
models. We find that our value of 8q is consistent with their perturbative determination 
for the 0(4) case, presented in [T2j. However, we do not confirm their conjecture that 



fit (a + dj) 


9 P 


Zp 


Jx\ z p) 


3+1 


0.605(2) 


1.25(1) 


0.345(2) 


2+2 


0.604(3) 


1.21(3) 


0.344(2) 


1+3 


0.557(3) 


1.33(2) 


0.3553(7) 



Table 10: Results for 9 P , z p and f x (z) using fits reported in Table 



the Cj's get smaller with increasing i. 

We also stress that, in addition to providing a better fit to the numerical data, the 
expression considered is a continuous function, needing no interpolation between the two 
x regions. This is particularly useful for the determination of the pseudo-critical line, 
since the interpolating form introduced in ^U] is unstable precisely in this region. As a 
result, our determination of z v is very precise in comparison to the previous estimates 
from the interpolated form and the perturbative equation of state. 
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